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ABSTRACT 

In a recent publication, the flexion aperture mass statistic was found to provide a ro- 
bust and effective method by which substructure in galaxy clusters might be mapped. 
Moreover, we suggested that the masses and mass proflle of structures might be con- 
strained using this method. In this paper, we apply the flexion aperture mass technique 
to HST ACS images of Abell 1689. We demonstrate that the flexion aperture mass 
statistic is sensitive to small-scale structures in the central region of the cluster. While 
the central potential is not constrained by our method, due largely to missing data 
in the central 0.5' of the cluster, we are able to place constraints on the masses and 
mass proflles of prominent substructures. We identify 4 separate mass peaks, and use 
the peak aperture mass signal and zero signal radius in each case to constrain the 
masses and mass proflles of these substructures. The three most massive peaks exhibit 
complex small-scale structure, and the masses indicated by the flexion aperture mass 
statistic suggest that these three peaks represent the dominant substructure compo- 
nent of the cluster (~ 7 x lO^^/i^^M©). Their complex structure indicates that the 
cluster - far from being relaxed - may have recently undergone a merger. The smaller, 
subsidiary peak is located coincident with a group of galaxies within the cluster, with 
mass ^ 1 X 10^^h~^MQ. These results are in excellent agreement with previous sub- 
structure studies of this cluster. 
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1 INTRODUCTION 

Accurately characterising the mass profiles of large struc- 
tures in the universe yields important constraints on our cos- 
mological model and theories of structure formation. While 
N-body simulations carried out under the assumption of the 
standard concordance cosmological model suggest that cold 
dark matter haloes are well fit by Navarro- Frenk- White mass 
density profiles (Navarro, Frenk & White, 1997), gravita- 
tional lensing observations of field galaxies have found to- 
tal density profiles that are consistently close to isothermal 
(Treu & Koopmans, 2002; Rusin, Kochanek & Keeton, 2003; 
Rusin & Kochanek, 2005; Koopmans et al., 2006; Gavazzi 
et al., 2007; Parker et al., 2007; Czoske et al., 2008; Dye et 
al., 2008; Tu et al., 2009), refiecting the impact of baryonic 
physics. The NFW model has only recently become favoured 
over the isothermal model in studies of clusters of galaxies, 
away from their inner regions (e.g. Carlberg et al., 1997; van 
der Marel et al., 2000; Athreya et al., 2002; Lin et al., 2004; 
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Hansen et al., 2005; Katgert, Biviano & Mazure, 2005; Lokas 
et al., 2006; Rines & Diaferio, 2006; Wojtak et al., 2007; Ok- 
abe et al., 2009). 

Moreover, whilst N-body simulations suggest an anti- 
correlation between mass and NFW concentration parame- 
ter in massive haloes (e.g. Navarro et al., 1997; Shaw et al., 
2006), and such a correlation has been seen in lensing stud- 
ies of field galaxies (Mandelbaum et al., 2008), some clusters 
of galaxies appear to have much higher concentrations than 
would be expected. One example of such a cluster is Abell 
1689, where previous strong and weak lensing studies have 
found best-fit mass profiles that diflFer from one another by 
up to a factor of 4 in mass and 3 in concentration parameter 
(see, for example, Peng et al., 2009 and references therein). 

Measurements of fiexion oflFer an independent method 
by which the density profiles of clusters of galaxies and their 
associated substructures might be characterised (Goldberg 
& Bacon, 2005; Bacon et al., 2006). Indeed, fiexion mea- 
surements - being a direct probe of the gradient of the lens- 
ing convergence - might well allow a clear discrimination 
between mass profiles. In a recent publication, Lasky and 
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Fluke (2009) have shown that the first flexion signal from 
an SIS profile differs from that of an NFW profile substan- 
tially at moderate angular separation between the source 
image and the lens, whilst the second flexion signal shows 
strong variation when the NFW concentration parameter is 
varied. 

In Leonard, King & Wilkins (2009), an aperture mass 
statistic for gravitational flexion was derived and shown to 
provide a robust method by which massive structures may 
be detected and mapped out. Moreover, it was shown to be 
more effective at detecting small-scale substructure than tra- 
ditional weak lensing methods. Further exploration of this 
method (Leonard & King, 2010) demonstrated, with sim- 
ulated data, that multiple aperture mass reconstructions, 
performed under differing aperture parameters, might allow 
one to determine the mass and mass profile of a structure 
or substructure in an entirely nonparametric way, even with 
relatively low signal to noise in the reconstructions. 

In this paper, we demonstrate this method using flex- 
ion data from HST ACS images of Abell 1689. These data 
were previously used in a galaxy-galaxy flexion study of the 
cluster, as well as to generate a parametric model of the sub- 
structures in the galaxy cluster (Leonard et al., 2007). In §(2] 
we briefly review the weak lensing and flexion formalism, the 
flexion aperture mass statistic FM^p , and describe a method 
for model discrimination using FM^p reconstructions. In §[3l 
we review previous lensing studies of Abell 1689, and out- 
line the data and flexion measurement methods. In § [4] we 
present the results of both our flexion aperture mass recon- 
structions, and place constraints on the locations, masses, 
and mass profiles of the substructures detected by the FM^p 
method. We discuss our results in the context of previous 
studies of the cluster in § O and conclude with a summary 
of our findings in § [B] 

Throughout the text, we assume a standard ACDM 
cosmology with Qm — 0.27, Ha — 0.73, and Hq — 
100/ikms"^Mpc"\ 



2 REVIEW OF FORMALISM 

In this section, we briefly review the formalism underlying 
traditional weak lensing and flexion studies, and present the 
aperture mass statistic for flexion. While this statistic does 
not directly provide a determination of the masses of struc- 
tures detected, it does provide a robust, non-parametric 
method for detection of mass concentrations (see, e.g. 
Schneider, 1996; Schneider et al., 1998; Leonard, King & 
Wilkins, 2009). 
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is the magnification matrix, which is related to the dimen- 
sionless projected surface mass density k of the lens and the 



complex shear 7, and the D operators Dijk 
given by (Bacon et al., 2006) 
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where T and Q are first and second flexion, respectively. 

While the convergence gives rise to an isotropic mag- 
nification of a background image, the shear and flexion op- 
erators yield anisotropic distortions in galaxy images. The 
complex shear 7 = 71 +^72 = 1716^*"^ induces a tangentially- 
aligned ellipticity in galaxy images, while the first flexion, 
T — T\ ^ iT2 — \T\e^'^ ^ is a direct probe of the gradient 
of the convergence, and gives rise to a skewness and a shift 
in the centroid of the brightness distribution of the galaxy 
image. Second flexion, Q — Q\^iQ2 — \Q\^^'^ ^ has m = 3 ro- 
tational symmetry, and gives rise to a bending or arciness in 
galaxy images. As described in Goldberg & Leonard (2007; 
see also Leonard et al., 2007; Okura et al., 2007, 2008), reli- 
able measurements of this latter distortion have not, to date, 
been obtained; therefore for the remainder of this paper, we 
focus on measurements of first flexion only. 

Each of these distortions may be directly related back 
to the convergence through a convolution of the form 
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where Me is the component of the shear (flexion) aligned 
tangentially (radially) with respect to the centre of the lens, 
and the convolution kernels are given by 
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Lasky & Fluke (2009) have presented analytic expressions 
for the convergence, shear and flexion for a number of 
commonly-used mass profiles. Relevant results for the sin- 
gular isothermal sphere and NFW models considered in this 
paper are given in Appendix [Al 



2.1 Shear and Flexion 

First, we begin by describing the origin of the shear and 
flexion terms in weak lensing. As described in Goldberg & 
Bacon (2005), in the case of a smoothly varying lens field, the 
lens equation describing the mapping between coordinates 
in the source plane, and coordinates in the lens plane, 
Oi, can be expressed as 

Pi — AijOj + —DijkOjOk , (1) 



2.2 The Flexion Aperture Mass Statistic FMap 

The flexion is a method by which the first flexion is related to 
the underlying lens convergence. In this method, apertures 
of radius R are laid down at various points over the field. 
Within each aperture, the convergence is convolved with a 
mass filter function w{r). This convolution is then related 
to the measured flexion convolved with an appropriate filter 
function Q(r); i.e. for an aperture located at angular position 
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Figure 1. Mass {left panel) and flexion {right panel) fllter func- 
tions as described in the text for i =3, 5, 7, and 10. 



xo, the flexion aperture mass statistic is defined as 
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where x is the co-ordinate with respect to the centre of the 
aperture and Te is the compontent of the flexion oriented 
radially with respect to the aperture centre. Note that the 
definition of FM^p given by equation [8] is formally identical 
to the definition of the shear aperture mass statistic (see 
Schneider 1996). The two methods are distinct only in the 
measurable quantity considered and the choice of filter func- 
tion Q, which is related to the mass filter function w (see 
equation [11] below) . 

The mass filter function is subject to the condition 
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in order to ensure that the aperture mass statistic is invari- 
ant under a transformation /^(x) + ^o, where kq is 
an additive constant. Further, the mass filter function, w, 
must go to zero smoothly at the aperture boundary. The 
filter functions w and Q can be related to each other by 
considering the relationship between the convergence k and 
the flexion (see Leonard, King & Wilkins, 2009) 



1 /""^ / / / 
Q{x) — / x'w{x') dx' 

X Jo 



w{x) 



-Q[x) 



dx 



(11) 



In this paper we choose our mass filter function w to be 
given by a family of polynomial functions (Schneider et al., 
1998) 
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is a normalisation constant, and I is an integer defining the 
polynomial order of the filter function. This family of mass 
filter functions yields the flexion filter functions 
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These filter functions are plotted in Figure[T]for various val- 
ues of I. 



As discussed in Schneider et al. (1998) and Leonard, 
King & Wilkins (2009), these filter functions are not de- 
signed to provide optimal signal to noise for any given mass 
profile; rather, the advantage of using such a family of filter 
functions is that they provide a robust way to detect real 
structures from lensing measurements. Moreover, higher- 
order polynomials provide more sharply peaked filter func- 
tions, which are sensitive to structures on smaller scales. 
Therefore, changing the aperture size and the polynomial 
order will provide information on a range of different phys- 
ical scales. 

Filter functions that better match particular shear pro- 
files have been suggested in the literature; for example 
Schirmer (2004) proposed a filter more appropriate to the 
NFW profile in the context of searching for galaxy clus- 
ters using weak lensing shear data. In addition, Bartelmann, 
King & Schneider (2001) considered the impact of lens el- 
lipticity on the measured shear aperture mass (Map) signal, 
concluding that for a moderate ellipticity there is only a 
few percent difference compared with a spherically symmet- 
ric lens. Similar corrections would be expeted for the FMap 
statistic, as it is formally identical to that for shear. 

One of the advantages of using aperture mass tech- 
niques is that the noise in the reconstruction is fairly 
straightforward to estimate, and the noise properties of 
the aperture mass maps are generally well-understood. The 
noise may be estimated by repeatedly rotating the individ- 
ual flexion vectors through a random angle, and repeating 
the measurement. The standard deviation of an ensemble of 
such randomisations gives a good estimate of the reconstruc- 
tion noise. The mean signal from these randomisations also 
gives an indication of the level to which systematic errors 
affect the measurement. Zero signal is expected; deviation 
from this indicates a systematic effect. For a more complete 
discussion of the noise properties of the aperture mass statis- 
tic, the reader is referred to Schneider (1996). 

Systematic errors may also be probed by considering 
the so-called B-mode signal, obtained by rotating the flex- 
ion measurements by 7r/2 radians and repeating the aperture 
filtering. Gravitational lensing naturally does not give rise 
to B-modes; therefore a B-mode signal in the resulting aper- 
ture mass maps is indicative of residual systematic errors in 
the measurement, e.g. arising from incomplete correction of 
galaxy ellipticities for the impact of telescope optics. 

2.3 Determining Masses and Profiles with FM^p 

In Leonard & King (2010), we considered the radially- 
averaged flexion aperture mass signal from isothermal and 
NFW profile lenses with varying concentration parameter, 
and investigated the evolution of the signal under changes to 
the aperture parameters (the radius R and polynomial order 
i). Specifically, measurements of the peak signal and zero- 
signal radius of the aperture mass signal for various aperture 
parameters allowed us to discriminate between mass models 
using simulated flexion data. In this paper, we apply this 
method to ACS images of Abell 1689. 

The method centres on comparing data to simulated 
template profiles, and cross-matching data and templates 
as the aperture parameters are varied. For each aperture 
mass reconstruction, mass peaks are identified in the signal 
to noise map. The aperture mass signal is then circularly- 
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averaged around the centre of the signal-to-noise peak, and 
the peak signal and zero-signal radius are computed. The 
template models are then searched, and those models show- 
ing a peak signal and a zero-signal radius that match the 
data, within the error bars, are identified. As the aperture 
parameters are varied, the set of template models compat- 
ible with all the data points is reduced, as a result of the 
fact that the peak signal and zero signal radii from different 
mass models vary differently under changes to the aperture 
parameters. 

The template profiles span a range of masses between 
10^^ - 2 X lO^^/i"^M0, and use both SIS and NFW models 
with concentration parameter sampled at integer values in 
the range c = 3 — 12. The template models are generated 
assuming a constant background source redshift of Zq = 1.0 
and a lens redshift of za = 0.2. 

These distances come into the normalisation of the flex- 
ion signals in different ways for the SIS and NFW density 
profiles. As described in Leonard & King (2010), the SIS 
model induces a flexion related to the source and lens sepa- 
rations through 

^sis oc ^ , (15) 

whilst the redshift dependence of the NFW flexion signal is 
given by 

J'nfw oc H'^{zd)^^^ , (16) 

where H{zd) = ^^o\/^m(1 + ^d)^ + ^a- By appropriately 
scaling the template profiles, it is straightforward to com- 
pare measurements with different source and lens redshift s. 

One might imagine that in assuming a constant back- 
ground source redshift, we might introduce a bias in our 
measurements, which may manifest as a tilt in the radial 
profile of the aperture mass statistic due to magnification 
effects in the neighbourhood of mass concentrations. How- 
ever, we ignore this effect for the purposes of this paper. 
Given that the correction factors described above are fairly 
small under changes to the source and lens redshift s, and 
given that the FM^p reconstructions themselves are rather 
noisy, it is unlikely that this bias will be a dominant source 
of error. 



3 ABELL 1689 

3.1 Gravitational Lensing Studies 

Abell 1689 is a massive cluster of galaxies located at a red- 
shift of ^ = 0.187 (Frye et al., 2007), which exhibits the 
largest known Einstein radius for a gravitational lens ob- 
served to date: Oe = 45^' for Zs = 1 (see, e.g., Tyson, Wenk 
k Valdes., 1990; Miralda-Escude & Babul, 1995; Broadhurst 
et al., 2005a, b). This cluster has been the subject of numer- 
ous X-ray and gravitational lensing studies, and the mass 
estimates derived from these two methods have often been 
found to disagree by a factor of 2 or more (see Peng et 
al., 2009, and references therein). The X-ray morphology 
of the cluster appears to be fairly regular, implying hydro- 
static equilibrium. However Andersson & Madejski (2004) 
find an asymmetric temperature distribution as well as a 



high-redshift substructure within the cluster, which suggests 
that Abell 1689 is currently undergoing a merger. 

The existence of substructures, contributing ~ 7% of 
the total mass within 250/i~^ kpc (~ 1-34^), has been con- 
firmed using gravitational lensing (e.g. Broadhurst et al., 
2005b; Diego et al., 2005; Halkola, Seitz, & Pannella, 2006; 
Zekser et al., 2006; Limousin et al., 2007; Saha et al., 2006), 
though these substructures may be extended structures pro- 
jected along the line of sight but external to the cluster, 
rather than merging clumps and groups within the cluster 
system. Indeed, Lokas et al. (2006) have studied the redshift 
distribution of galaxies in the cluster and its surroundings 
and concluded that there are structures along the line of 
sight that will affect lensing, but not dynamical, mass es- 
timates. Such structures could help to explain the discrep- 
ancy in mass estimates between those derived through lens- 
ing measurements and those derived from X-ray studies, the 
latter of which tend to be systematically lower. 

Even amongst gravitational lensing studies themselves, 
however, discrepancies are seen. Strong lensing mass esti- 
mates appear to be higher than those obtained using weak 
lensing alone, with values ranging from M200 = 0.85 x 
10^^ Mo (King, Clowe, k Schneider, 2002) using weak lens- 
ing to M200 = 3.05 X 10^^ M© (Halkola et al., 2006) using 
strong lensing. Furthermore, gravitational lensing studies of 
Abell 1689 and other clusters often show higher concentra- 
tion parameters than predicted by the mass-concentration 
relation seen in N-body simulations (e.g. Shaw et al., 2006), 
and show a wide spread in values across the lensing mea- 
surements. 

It has been argued that using a prolate halo model 
aligned along the line of sight, rather than a spherically- 
symmetric model, can help bring the measured concentra- 
tion parameter more in line with that predicted by simu- 
lations (Oguri et al., 2005; Corless & King, 2008; Corless, 
King, & Clowe, 2009). Indeed, Shaw et al. (2006) find that 
higher mass haloes tend to be prolate, rather than spheri- 
cal. Moreover, Morandi, Pedersen & Limousin (2010) have 
found that taking proper account of the triaxiality of DM 
haloes can help resolve the discrepancy between lensing and 
X-ray mass estimates of clusters. However, such models have 
thus far not reduced the measured concentration parameter 
appreciably in studies of Abell 1689. 

In addition, there is little agreement regarding the im- 
portance of substructure in the mass models used to fit the 
lensing data. For example, Limousin et al. (2007) found it 
necessary to include in their strong lensing model a sec- 
ondary mass peak, offset from the central core, in order to 
best fit the strong lensing features. On the other hand, oth- 
ers (e.g. Broadhurst et al. 2005a) have found a good fit with 
a single NFW halo. 

Most recently, Coe et al. (2010) have presented a com- 
bined strong and weak lensing analysis of the cluster, provid- 
ing the highest resolution mass reconstruction of this clus- 
ter to date, and using the same ACS images considered here, 
supplemented with ground-based imaging. Their reconstruc- 
tion shows complex substructure, and they find a best fit 
NFW profile with M200 = 1-8 x lO^^h'^MQ and c = 9.2 
from combined strong and weak lensing constraints. 

It is clear that Abell 1689 is a cluster with a complicated 
structure, and that there is still much work to be done to 
understand and accurately map the structure of the cluster. 
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As flexion directly probes the gradient of the convergence, 
it is sensitive to structure on smaller scales than traditional 
weak lensing studies and can probe a greater distance from 
the centre of the cluster than strong lensing studies are able 
to. For this reason, flexion studies might be of great impor- 
tance in the study of Abell 1689 and similar clusters. 

3.2 HST ACS Data and Analysis Pipeline 

The data consist of 20 HST ACS images of Abeh 1689, taken 
using the WFC. Each image has been reduced through the 
HST ACS pipeline and contamination from cosmic rays re- 
moved. These 2300-2400 second exposures cover a square 
field of view of 3.4^ on a side, have an angular resolution 
of 0.05''pixer\ and are described in full detail in Broad- 
hurst et al. (2005a). Of these images, 4 were taken using the 
F475W filter (G-band), 4 using the F625W filter (R-band), 
5 using the F775W filter (I-band) and 7 using the F850LP 
filter (Z-band). The data analysis pipeline used in this pa- 
per is described in full in Leonard et al. (2007), and the key 
points are summarised below. 

3.3 Catalogue Generation 

For the purposes of source detection, stacking the images 
within each colour band using the SWarp software package 
(Bertin et al., 2002) was found to be beneficial. The software 
first estimates and subtracts off the background in each indi- 
vidual frame. Then the images are resampled and projected 
on to the output frame using a tangential (gnomonic) pro- 
jection. We opted to use the median stacked image, rather 
than the mean, in order to avoid any biasing of the output 
image by spurious hot spots and bad pixels. A single co- 
added image was created for each of the G-, R- and I-bands. 
In the Z-band (as discussed in Leonard et al., 2007) it was 
found that very little benefit was gained in source detection 
from stacking 7 images as opposed to 3 or 4; therefore, two 
co-added images were created. 

Source extraction was carried out in two stages using 
SExtractor (Bertin & Arnouts, 1996), using a modification 
of the "hot" and "cold" strategy of Rix et al. (2004). During 
the first pass, a catalogue of known foreground objects in- 
cluding stars and cluster members either listed in Due et al. 
(2002) or listed as such in the NASA/IPAC Extragalactic 
Database (NED) was provided, and only those objects de- 
tected. These sources were then masked out by setting their 
pixel values equal to random noise, making use of the vari- 
ance maps produced by SExtractor, and thus simulating an 
emptier field. SExtractor was then run on the masked im- 
ages using a lower detection threshold, in order to detect 
the smaller, fainter background galaxies of interest. A final 
catalogue of objects was generated including all those ob- 
jects that were detected in at least two co-added images, 
including either the 1- or R-band image. 

The stacking process described above clearly will affect 
the PSF of the resulting median-stacked image in a complex 
way; each frame is offset and rotated slightly from the oth- 
ers within the stack, therefore the resulting PSF will be im- 
possible to model simply. However, Leonard et al. found the 
effect of the telescope point spread function to be significant 
only for shear measurements, but not for flexion measure- 
ments, provided the source is well-resolved. This is because 



the induced flexion resulting from the PSF is related to the 
image size by (see Leonard et al 2007 for full details of this 
calculation) : 

4 

^induced ~ ^PSF ^ — ^ , (17) 

^source ~r ^^psF 

where J^psf is the flexion measured in the PSF (typically 
of order 10~^ — 10~^/'^), apsF is the semi-major axis of 
the PSF (~ 0.1'') and ttsource is the semi-major axis of the 
source. Therefore flexion measurements were carried out on 
the stacked images, which offer improved signal to noise, 
and all sources with asource < apsF were rejected from the 
analysis. Flexion measurements were carried out using the 
HOLICs technique (Okuraet al., 2007; Goldberg & Leonard, 
2007). 

In addition to PSF effects, images from the ACS WFC 
suffer from a geometric image distortion resulting from the 
off-axis location of the camera on the Hubble Space Tele- 
scope. This distortion has been modelled by Meurer et al. 
(2003), and is weh fit by a 4*^ order polynomial. This fit 
was used to evaluate the shear and flexion induced by the 
geometric distortion as a function of location on the chip. 
Leonard et al. (2007) found that here, again, the effect on 
the flexion was negligible (^ 10""^/''), and note also that 
the geometric distortion is greatly reduced in the stacked 
frames as a result of the reprojection of the images during 
the stacking process. 

Sources were rejected at various points in the analysis 
pipeline, in order to reduce contamination from erroneous 
noisy measurements and to ensure that only well-resolved 
sources were used. Full details of the selection criteria can 
be found in Leonard et al. (2007). The resulting background 
source count is ^ 75arcmin~^. As in Leonard et al. (2007), 
we assume a median source redshift of Zs — 0.9, in accord 
with the measurements of Limousin et al. (2007) 

3.4 Aperture Mass Method 

The aperture mass statistic was computed using the poly- 
nomial filter described above, and using 1000 x 1000 evenly 
spaced apertures across the field. Due to the high resolu- 
tion of the images and the resulting high background source 
count, including a large number of randomisations in the 
noise calculation proved rather time-consuming. Therefore, 
only 100 randomisations were carried out for each aperture 
mass reconstruction. This was found to be sufficient, how- 
ever; increasing the number of randomisations for a given 
reconstruction was not found to significantly alter the error 
estimates. 

FMap reconstructions were carried out for aperture radii 
R = [45'', 60'', 75", 90"] and filter functions with polynomial 
order I = [3,5,7,10]. This set of aperture parameters will 
allow us to map the structure in the cluster on a range of 
scales, and provide a consistency check to help eliminate 
peaks arising purely out of noise. Moreover, having 16 dif- 
ferent reconstructions should allow to discriminate between 
mass models for any structures detected; Leonard Sz King 
(2010) demonstrated that, particularly when probing low- 
mass structures, a larger number of reconstructions is helpful 
in discriminating between competing mass profiles. For each 
reconstruction, a B-mode map was also generated in order 
to assess the level of systematic errors within the maps. 
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In order to more quantitatively characterise the struc- 
tures detected, an aggregate signal to noise map for FM^p, 
consisting of the average signal to noise of all 16 reconstruc- 
tions in each case, is generated. Such a map is useful for 
identifying persistent structures across the 16 reconstruc- 
tions, and is used solely for this purpose. As discussed in 
detail in Schirmer (2004) in the context of shear, significant 
Map peaks that are in fact due to noise disappear at filter 
scales only slightly larger or smaller than the scale where 
they are detected. In contrast, real peaks are detected as 
significant over a much broader range of filter scales. A sim- 
ilar argument will apply for the FMap statistic. 

Using the CLFIND algorithm developed by Williams, 
de Geus & Blitz (1994), we identify all significant peaks 
within the map down to a limiting signal to noise ratio. 
This software allows us to identify peaks even where they 
occur as a sub-peak in a larger structure. This is useful where 
structures have been broken up into smaller, blended sub- 
structures in the reconstructions. A location for each peak 
is then identified by computing the signal-to-noise weighted 
centroid. 

Similarly, for each individual reconstruction, peaks are 
identified and cross-matched with the aggregate map. A ra- 
dial profile for each is computed by averaging the signal in 
radial annuli about the signal-to- noise weighted centroid. 
Note that the centroid used for a given peak is not fixed; 
i.e. the centroids for each peak are independently calculated 
in each FMap reconstruction, and the radial average is com- 
puted with respect to that centroid, rather than the centroid 
computed in the aggregate map. In addition, a B-mode ra- 
dial profile is computed for each peak. Radial profiles ex- 
hibiting a B-mode signal not consistent with zero, within 
the error bars, are rejected from further analysis. 

3.5 Model-selection and mass determination 

As described in Leonard & King (2010), the radial profiles 
of the FMap signal can be used to discriminate between 
competing mass profiles. For each radial profile computed, 
a polynomial function is used to fit the signal, and the first 
real, positive root taken to be the zero-signal radius Rq. The 
radial profile of the associated noise map is also computed, 
and this yields error bars on each point in the radial profile 
for a given peak. 

The radial profiles can then be compared with the tem- 
plate models, and compatible regions of parameter space 
isolated for each peak. The catalogue is first searched for 
models whose peak signal falls within the error bar associ- 
ated with the measured peak of the radial profile. Models in 
this subset are then excluded if their zero-signal radius falls 
outside the range allowed by the measurement given the er- 
ror bars. The allowed range is taken to be all those values 
of xo for which the signal is compatible with zero given the 
1 — a error bars (see Figure[2|). In a number of cases, par- 
ticularly those in which the peak has substantial associated 
substructure, this only allows us to place a lower limit on 
the zero-signal radius. 

This is a very conservative constraint, allowing a large 
range of Rq values to be considered in each case. However, 
given the low signal to noise seen in our maps (see §[4lbelow), 
it seemed prudent to allow a larger range of Rq, so as not 
to present overly optimistic results. 



Peak 1, I = 5, R = 60 Peak 2, I = 5, R = 60 




Figure 2. Sample radial profiles for each of four peaks detected in 
the FMap reconstruction with £ = 5 and R = 60" (see §|4] below) 
showing the FM^p signal and associated error bars, the B-mode 
signal, the polynomial fit used to determine Rq (solid curve), and 
the range of Rq considered to be compatible with the data in the 
subsequent analysis (indicated by vertical dashed lines). 



Figure[2] shows four example radial profiles to illustrate 
the method described above. In this figure, we plot the 
radially- averaged FM^p signal for each of four peaks identi- 
fied in the FM^p reconstruction with ^ = 5 and R = 60^^ and 
the associated error bars. We overlay the polynomial fit used 
to determine the value of Rq as a solid line. Dashed vertical 
lines indicate the range of Ro considered to be compatible 
with the data; the range of mpeak considered is that covered 
by the xo = error bar. We also plot the radial profile of 
the B-mode, to demonstrate the criteria for inclusion in the 
subsequent analysis. 

The left panels show examples of peaks for which we 
were only able to place a lower limit on the range of Rq 
compatible with the data, whilst the right panels highlight 
examples for which a smaller range of Ro is considered. In 
addition, the figure indicates which peaks are included or 
rejected based on their B-modes. The peaks in the right 
column clearly do not exhibit appreciable B-modes and are 
included. The profile shown in the top left panel of the figure 
has B-modes that are clearly not consistent with zero within 
the error bars, and is therefore rejected. The peak in the 
bottom left panel again shows an appreciable B-mode signal; 
however the signal is compatible with zero within the error 
bars for all values of xo, and therefore this peak is included in 
the subsequent analysis. Such borderline cases are relatively 
rare though, occurring for approximately 1 out of every 6-7 
profiles considered. 

For each model included for a given value of i and i?, 
we compute a fractional deviation from the measured signal. 
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given by: 



2 _ (^peak — '/^p, model) ^ (^0 — i?0, model) 



peak 



(18) 



After iterating over the available reconstructions, a weight 
W can be assigned to each model 



where 



(A) ' 



(19) 



(20) 



is fractional deviation for the model averaged over all realisa- 
tions for which it was identified. In effect, (A) gives us a sort 
of mean statistic in parameter space for a given model 
as compared with the data. More accurately, (A) gives the 
mean distance (in parameter space normalised by the mea- 
sured FMap parameters) of a given model from the data. 
The statistic ensures that greater weight is given to mod- 
els for which A^det, the number of realisations for which the 
model was identified as compatible with the data, is higher 
and for which the mean deviation between the model and 
the data is smaller. 

Note that we treat the E-mode data for those profiles in- 
cluded in subsequent analysis entirely uncritically. In other 
words, no down- weighting is applied to those cases which 
have marginal B-modes as compared to those with none. 
This will introduce errors in the subsequent analysis; how- 
ever we do not expect this to bias our results significantly, 
as such borderline cases are relatively rare. 

This method allows us to identify models that are com- 
patible not only with the peak signal and zero signal radius 
in a given realisation, but also with the evolution of these 
parameters as the aperture parameters are varied. So, while 
a simple x^-fitting of a class of model (such as the NFW 
model) to the data will tell us how well a given model fits 
a given radial profile at fixed I and i?, the W statistic tells 
us how well a given model fits the ensemble of FM^p radial 
profiles. In this way, W acts as an unnormalised pseudo- 
likelihood distribution, yielding higher values for models 
which better fit the data, and are therefore more likely rep- 
resentative of the true matter density profile. 

Whilst this is a non-traditional statistic, and - being 
unnormalised - may be seen as somewhat difficult to inter- 
pret, there are several advantages to this approach. Firstly, 
a full Bayesian treatment of the flexion aperture mass statis- 
tic requires both a knowledge of the intrinsic distribution of 
flexion in unlensed galaxies, and the response of this dis- 
tribution to aperture mass filtering on different scales. To 
date, the intrinsic flexion distribution in field galaxies not 
been measured or characterised, but a simple gaussian ap- 
proximation is unlikely to suffice (Rowe et al. 2007), as is 
usually taken to be the case for shear. 

In addition, while Bayesian methods yield straightfor- 
ward error estimates on the best-fit parameters (which are 
not obtained using our method), they do not allow to com- 
pare between different families of models. In other words, 
a maximum-likelihood estimator can compare NFW mod- 
els with different masses and concentration parameters, but 
cannot directly compare that class of models to the class 
of SIS models simultaneously. Our method allows for such a 



comparison - by defining the W statistic in exactly the same 
way for the SIS and NFW models, we are able to quantita- 
tively compare the goodness-of-fit of a given SIS model with 
an NFW. 

This is not to imply that standard statistical estima- 
tors are not important. In fact, in order to fully quantify 
the goodness-of-fit of a given model, parametric modelling 
will be essential. What we present here is simply a fast and 
effective method to identify substructure, to discriminate be- 
tween competing families of models, and to eliminate areas 
of parameter space incompatible with the data efficiently. 



4 RESULTS 

4.1 FMap reconstructions and peak identification 

The signal to noise contours of the 16 FM^p reconstruc- 
tions are plotted in FigurefS] overlaid on an image of the 
cluster. Here again, the contours are plotted at levels S = 
[1.0, 1.5, 2.0, 2.5]. The peak signal to noise value for each of 
the structures detected vary with each reconstruction, gen- 
erally being reduced slightly as the aperture size is increased 
at constant /, but is consistently found to be above <Speak 2 
for the larger, dominant structures. 

Several interesting features arise in the reconstructions. 
Firstly, as expected, the reconstructions using smaller aper- 
tures reveal structures on finer scales. Indeed, as expected, 
when the filter polynomial order is increased for a fixed 
aperture size, a similar effect is seen. For example, re- 
markable similarity is seen in the signal to noise plots for 
^ = 5, R = 45'' and ^ = 10, R = 60'' . This is an important 
consistency check of the method. The FM^p reconstructions 
are not independent as the filter functions themselves are not 
orthogonal. Therefore, any large inconsistencies between the 
maps, particularly those whose filter functions peak at the 
same angular scale, are indicative of persistent systematic 
or computational errors. 

It is clear that there are several prominent, persistent 
features across the 16 maps. These structures can be seen 
to break up into smaller substructures in the signal to noise 
maps as the aperture size is reduced or the polynomial order 
increased. In order to determine the most significant persis- 
tent structures in the FMap reconstructions, we consider the 
aggregate map, plotted in FigureH] 9 peaks are identified in 
total; for each, the signal-to-noise weighted centroid is high- 
lighted in the figure, and the peaks are labelled in order of 
decreasing signal to noise in the aggregate map. 

We now consider the peak detections in the individual 
FMap reconstructions. The signal to noise maps for each 
reconstruction are cross-matched with the average map, and 
where a peak has been broken down into smaller sub-peaks, 
the sub-peaks are identified by considering all peaks in the 
individual map that are located within the <S = 1 contour of 
the peak in the average map. 

In addition, a signal-to-noise weighted centroid is com- 
puted for the collection of sub-peaks corresponding to the 
primary peak in the average map. Figure[5] shows the peaks 
detected in each FMap realisation, appropriately labelled at 
the centroid of each peak or sub-peak. Overlaid on each map 
are contours showing the absolute value of the B-mode map 
for each reconstruction, normalised by the noise in that re- 
construction. Contours are plotted for |FMap g| = [Icr, 2a]. 
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Figure 3. The figure above shows fiexion aperture mass reconstructions for Abell 1689 using a polynomial fiexion filter function with 
£ = 3 (top row), ^ = 5 (second row), i = 7 (third row) and £ = 10 (bottom row), and an aperture radius of R = Ah" (first column), 
R = 60" (second column), R = 75" (third column) and R = 90" (fourth column). 



The la dispersion in the locations of each of the peaks 
are overplotted in FigureSl and provide an indication of the 
extent to which the locations of the peaks are offset from 
one another in the 16 FMap reconstructions. Offsets in the 
peaks can arise due to noise in the reconstruction, or as a 



result of a peak breaking up in to smaller sub-peaks; a small 
dispersion in the peak locations is indicative of a robust 
detection of the position a given peak. 

Of the 9 peaks identified, four (peaks 5, 7, 8, and 9) are 
located at the edge of the field. The B-mode maps appear to 
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Figure 5. The figure above shows the detected peaks in FMap reconstructions of Abell 1689 using a polynomial fiexion filter function 
with i = 3 (top row), i = 5 (second row), i = 7 (third row) and i = 10 (bottom row), and an aperture radius of = 45" (first 
column), R = 60" (second column), R = 15" (third column) and R = ^0" (fourth column). Peaks and their subpeaks are labelled in 
each reconstruction at the location of the signal-to-noise weighted centroid of the peak. Dashed contours show the absolute value of the 
B-mode signal, plotted at levels of Icr and 2cr. 
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Figure 4. The average signal to noise map of the 16 FMa.p re- 
constructions, showing all peaks identified with a signal to noise 
greater than 1. Overplotted aon the map are the locations of each 
of the 9 peaks identified, labels for each peak, and the la disper- 
sion in the location of the peaks across the 16 reconstructions 

show some structure around the edges of the field, thus we 
cannot consider these to be true detections. Even if they are 
indeed real features, due to their locations it is impractical 
to attempt to compare their radial profiles with those of 
the template models, given that they are located so close to 
the edge of the field; a complete circular average around the 
centre of each peak is impossible to obtain given the data. 

In addition, significant B-modes were detected in all 
FMap reconstructions of peak 4, as is evident from FigureO 
This peak can therefore also be considered to arise out of 
noise, and is excluded from further analysis. Therefore, we 
restrict our study to peaks [1,3], 2, and 6. Each peak is 
considered separately below. Note that peaks 1 and 3 are 
considered together, as it is likely that they form part of 
the same structure. Indeed, in many of the large R, low / 
reconstructions, peak 3 is undetected. 

4.2 Constraints on mass profiles 

4.2.1 Peaks 1 & 3 

Peak 1 is detected in all 16 FM^p reconstructions, whilst 
peak 3 appears in only 9 reconstructions, forming part of 
peak 1 for the remainder. It therefore seems likely that peak 
3 is a slightly lower-mass substructure located in the neigh- 
bourhood of peak 1 that becomes blended with peak 1 as 
the aperture radius is increased. Both peak 1 and peak 3 
are seen to break up into smaller parts in the high low 
R aperture mass reconstructions, which indicates a complex 
structure that is unlikely to be well-modelled by a spher- 
ical halo. Nonetheless, Leonard & King (2010) found that 
the peak signal and zero signal radius of the flexion aper- 



ture mass statistic could be used to estimate the masses of 
haloes even when they possess a moderate ellipticity. 

In the instances where peak 3 is clearly resolved, it does 
not show prominent B-modes within the error bars of the 
FMap reconstruction in each case. This implies that it is, in- 
deed, a real mass concentration within the cluster. Peak 1 is 
found to show very prominent B-modes in the 7 reconstruc- 
tions in which peak 3 is not resolved, which suggests that 
in these reconstructions, the signal from peak 3 might be a 
strong contaminant to the signal. Peak 1 shows prominent 
B-modes in a further 3 reconstructions, all corresponding 
to low-^, large- fllters. These B-modes most likely result 
from edge effects, as peak 1 is located close to the image 
boundary. 

Excluding those proflles with substantial B-modes, we 
are left with 6 reconstructions of peak 1 and 9 of peak 3, and 
we can attempt to constrain the mass distribution associated 
with these peaks from the remaining reconstructions. 

We first consider the W distribution of peak 1. Figure[6] 
shows the measured W distribution for the set of NEW tem- 
plate models (top panel) and SIS models (bottom panel) . In 
the NEW case, contours are plotted in M200 — c space at lev- 
els equal to [0.5, 0.7, 0.9] x the maximum W found amongst 
afl (both NEW and SIS) models. Dotted and dashed con- 
tours indicate A^det for each model. In the SIS case, the solid 
curve shows the W distribution as a function of mass, and 
the dashed curve shows A^det- Dotted horizontal curves show 
the W values corresponding to the contours in the NEW 
plot. 

The best-fit NEW model to the data is one with M200 = 
5 X IO^^H-^Mq and c = 11 (Wmax = 9.38, A^det = 6), 
and the best-fit SIS has M200 = 2 x lO^'^h'^MQ (Wmax 
8.82, A^det = 6). Thus, the data seem to favour the NEW 
model, and Eigure[6] shows a rather narrow range of masses 
favoured by the data. 

It is important to note that the range of allowed masses 
is similar in both the NEW and SIS cases. This is an impor- 
tant consistency check, as one would not expect the signal 
from low mass SIS profiles to be comparable to that from 
high mass NEW profiles, and vice versa. Despite the pres- 
ence of significant structure within each peak at small aper- 
ture radii and high filter polynomial order, and the large 
error bars in the measurements of rripeak and Ro, we are 
able to eliminate large regions of parameter space as incom- 
patible with the data, and isolate a small range of best-fit 
models for the data. 

We now consider peak 3, and the results are shown 
in Eigure[7] below. The best-fit NEW model to the data is 
one with M200 = 2 x 10^^h~^MQ and c = 10 (Wmax 
19.79, A^det = 9), and the best-fit SIS has M200 = 2 x 
IO^^H-^Mq (Wmax = 12.13, A^det = 7). Hcrc, again, the 
data seem to favour the NEW model, and Eigure[7] shows a 
narrow range of allowed masses. 

As expected, the peak 3 data favour a lower-mass halo 
than those of peak 1, and imply an overall mass for the sub- 
structure system of M200 5 — 7 x lO^^/i^^M©, consisting 
of two prominent sub-clumps. It is clear that the use of a 
non-circular model and the presence of wider-field data are 
needed to constrain the mass distribution more accurately. 
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Figure 6. The pseudo-likelihood distributions for NFW (top 
panel) and SIS (bottom panel) models for peak 1 generated as 
described in the text. In the top panel, filled contours indicate 
W at contour levels equal to [0.1,0.3,0.5,0.7,0.9] x Wmax, with 
the model corresponding to the peak W in the NFW distribution 
denoted by a white asterisk. The solid red curve shows W = 0. 
Corresponding levels are shown in the bottom panel as dotted hor- 
izontal lines. The dotted and dashed contours in the top panel in- 
dicate A/'detj represented in the bottom panel as the dashed curve. 



Figure 7. The pseudo-likelihood distributions for NFW (top 
panel) and SIS (bottom panel) models for peak 3. In the top 
panel, filled contours indicate W at contour levels equal to 
[0.1,0.3,0.5,0.7,0.9] X Wmax, with the model corresponding to 
the peak W in the NFW distribution denoted by a white aster- 
isk. The solid red curve shows W = 0. Corresponding levels are 
shown in the bottom panel as dotted horizontal lines. The dotted 
and dashed contours in the top panel indicate A^det? represented 
in the bottom panel as the dashed curve. 



4.2.2 Peak 2 

As before, but to a lesser extent, in the case of peak 2 we 
find a rather complicated structure, with the peak breaking 
down into 2-3 sub-peaks in four of the FM^p reconstructions. 
However, a B-mode signal appears in only one reconstruc- 
tion (^ = 3, R = 90"). Excluding this reconstruction, and 
treating peak 2 as a single halo, we compute the W pseudo- 
probability distribution as before. The results are shown in 
Figure[8] 



For this peak, we find a best-fit NFW halo with M200 = 
2 X lO^^/i-^M© and c 12, with Wmax 13.30 and A^det 
11. Consideration of the A^det contours in this case suggests 
that a lower concentration parameter (perhaps as low as 
8) may be compatible for the same M200. The best fit SIS 
model yields the same M200 and A^det, with Wmax = 11.25. 
Note that we see a tail in the NFW probability distribution 
towards high-mass, low-c models, which does not appear in 
the SIS distribution. This indicates that substructure may 
play a role in overestimating the zero-signal contour, thus 
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Figure 8. The pseudo-likelihood distributions for NFW (top 
panel) and SIS (bottom panel) models for peak 2. In the top 
panel, filled contours indicate W at contour levels equal to 
[0.1,0.3,0.5,0.7,0.9] X Wmax, with the model corresponding to 
the peak W in the NFW distribution denoted by a white aster- 
isk. The solid red curve shows W = 0. Corresponding levels are 
shown in the bottom panel as dotted horizontal lines. The dotted 
and dashed contours in the top panel indicate A^det: represented 
in the bottom panel as the dashed curve. 

giving rise to the bias. It is important to note, however, that 
this does not affect our overall mass estimate, even when the 
constraint on Rq is removed. 

4.2.3 Peak 6 

This peak is the lowest signal to noise peak considered here, 
and does not appear in the ^ = 3, = 90^^ reconstruction, 
indicating it is likely to be a relatively low-mass structure 
with a small physical size (or high concentration parameter) . 




M2oo(h ^M^^J 

Figure 9. The pseudo-likelihood distributions for NFW (top 
panel) and SIS (bottom panel) models for peak 6. In the top 
panel, filled contours indicate W at contour levels equal to 
[0.1,0.3,0.5,0.7,0.9] X Wmax, with the model corresponding to 
the peak W in the NFW distribution denoted by a white aster- 
isk. The solid red curve shows W = 0. Corresponding levels are 
shown in the bottom panel as dotted horizontal lines. The dotted 
and dashed contours in the top panel indicate A^det? represented 
in the bottom panel as the dashed curve. 

In contrast to the other peaks, peak 6 shows no complicated 
structure, and is likely to be well-fit by a circular model. 
Furthermore, it shows no B-mode signal in any of the 15 
FMap reconstructions in which it is detected. The pseudo- 
probability distribution for this peak is shown in Figure[9] 

We find that this peak is best fit by an SIS model with 
M200 = 1 X lO^^/i-^Mo (Wmax = 15.62, Ndet = 13. The 
best-fit NFW model has M200 = 9 x 10^^ H'^Mq, c = 12, 
with only slightly lower Wmax = 14.66 and A^det = 12. 
Thus, a highly concentrated halo with mass ^ 
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is clearly implied with a fairly narrow spread of allowed 
masses, and a clear preference to high concentration param- 
eters (perhaps in excess of the maximum c— 12 sampled by 
our template models). 

4.3 Comparison of best-fit models to data 

In order to demonstrate how well our 'best-fit' models, as 
identified by the W statistic, actually fit the data, it is help- 
ful to compare directly the measured radial profile with the 
best-fit SIS and NFW model for each peak. 

As a representative sample of our results, we plot in 
FigurefTOl the radial profiles for all peaks for I — [3, 5, 7, 10] 
and R — 45", with the exception of peak 2, for which we 
use R = 75^^ This is because this peak breaks down into 
smaller sub-peaks in the R = 45" reconstructions. In these 
cases, the radial profile strongly under-estimates the peak 
signal and central regions of the profile, as a result of the 
fact that the centroid used does not correspond to the peak 
signal-to- noise for any of the sub-clumps. As a result, these 
reconstructions have few (if any) models compatible with 
their values of mpeak and i?o, and are not included in the 
calculation of the W statistic. It is for this reason that A^det is 
substantially lower than the total number of reconstructions 
of this peak. 

For clarity, Table[T] below details the properties of the 
4 peaks discussed above in each of the 16 FM^p reconstruc- 
tions. 

It is clear from Figure[TO] that the best-fit models not 
only fit individual radial profiles well, but also track the evo- 
lution of the profile under variation of ^ or consistently, in- 
dicating that the method is very effective at model-selection. 



5 DISCUSSION 

It is important at this point to compare our results with 
other studies of this cluster, both to compare the perfor- 
mance of the flexion aperture mass statistic on small scales 
with other lensing techniques, and to ensure that prominent 
structures identified in previous works are correctly identi- 
fied by our method. 

Whilst there is a wealth of published work on Abell 
1689, we restrict our comparison to four works: 

• Leonard et al. (2007), included to ensure consistency in 
results obtained using two distinct flexion techniques on the 
same data. 

• Limousin et al. (2007) and Coe et al. (2010), both of 
which aim to quantify the substructure content of the cluster 
using strong and weak lensing measurements, and both of 
which make use of the same ACS dataset considered here, 
supplemented with wide field imaging for their weak lensing 
measurements. 

• Riemer-S0renson et al. (2009), who consider the sub- 
structures seen in strong and weak lensing maps of the clus- 
ter from Limousin et al. (2007), and compare these to high- 
resolution X-ray temperature maps from Chandra. 

As in Leonard et al. (2007), we find that the central 
mass density concentration is undetected in all of our flex- 
ion reconstructions. As discussed in the earlier work, this re- 
sults from the fact that our masking scheme strongly favours 



detection of background sources in the periphery of the im- 
age over those in the centre. The lack of data in this region 
significantly hinders any attempts to constrain the central 
density profile. 

However, the peaks 2 and 6 shown in this paper are in- 
deed detected in Leonard et al. (2007) - although slightly 
offset between the two maps - and the elongation in the 
Leonard et al. (2007) convergence map in the direction of 
the peak [1, 3] suggests a positive signal may have been seen 
in the earlier work. However, the parametric reconstruction 
carried out in that paper assumed that the flexion signal 
originated from known cluster members, and did not allow 
for the possibility of a significant dark clump, which the peak 
[1, 3] system appears to be. Similarly, the small offset in the 
detected positions of peaks 2 and 6 are unsurprising. There- 
fore, our results are entirely compatible with the earlier work 
and, moreover, offer the advantage that no assumptions were 
made regarding the locations of substructure peaks within 
the cluster. 

Peak 2 also appears as a prominent feature in the 
parametric strong lensing reconstruction carried out by 
Limousin et al. (2007). Modelling this clump as an isother- 
mal sphere, they find a best fit mass of M200 = (1.25 ± 
0.3)timesl0^^h~^ Mq within a radius of 1^56, in good agree- 
ment with our constraints. As their model only constrains 
the central ^ 2.6^ of the cluster, their mass model does not 
include our peak [1,3] system or peak 6. However, given 
that both these systems lie well outside the critical region 
of the cluster, we do not expect a strong-lensing model to 
be particularly sensitive to these structures. 

In Riemer-S0renson et al., the lensing reconstruction 
performed by Limousin et al. (2007) is compared to a high- 
resolution X-ray map obtained using Chandra. The temper- 
ature profile of the cluster exhibits a substructure to the 
northeast of the centre of the cluster, in the same direction 
as, but rather offset from, our peak 2. They argue that the 
X-ray and lensing substructures seen are coincident, with 
the gas mass amounting ^ (9.7 ± 0.3) x 10^^ M©. 

No other significant substructure is seen in the X-ray 
map of the region covered by the ACS imaging. However, 
this is not entirely surprising given the coarser angular res- 
olution of the Chandra images and the relative insensitivity 
of the X-ray temperature map to physically small, low-mass 
structures. 

Coe et al. (2010) present an excellent, high resolu- 
tion, non-parametric reconstruction of the cluster using their 
LensPerfect software. Their reconstruction (c.f. their Fig- 
ures) shows a clear detection of our Peak 2, as well as 
smaller sub-clumps in the direction of (but offset from) our 
peak 6. They also appear to show a marginal detection of 
our peak 3, but no detection of peak 1. It is likely that this 
arises due to a lack of strong lensing data in this region of the 
image, which is substantially offset from the critical, strong 
lensing regime of the cluster. 



6 SUMMARY 

In this paper, we have investigated in detail the complex 
internal structure of the galaxy cluster Abell 1689, using 
measurements of the flexion aperture mass statistic. While 
the shear aperture mass statistic is seen to be strongly con- 
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Table 1. This table summarises the properties of each peak in each reconstruction in which it was detected with a signal to noise 
«5peak > 1- The table lists the signal to noise value at the centroid used to compute the radial profile for each peak reconstruction, 
highlights those reconstructions exhibiting significant B-modes, and indicates whether the best-fit model exhibits values of mpeak and 
i?o compatible with the data. The right-most column lists possible reasons that a peak with no noticeable B-mode signal might not have 
been included in the final AT^et for the best-fit model. 





Figure 10. The figure shows a representative sample of radial profiles for the four detected peaks. Each row represents a single peak, i 
increases left to right, and R is fixed in each row. the best-fit models are overplotted as dashed red lines. 



taminated by systematic errors, the flexion aperture mass 
statistic appears to hold a great deal of promise for substruc- 
ture studies within clusters of galaxies. In the central 3.5^ of 
the cluster, we identify three significant structures, two of 
which break down into smaller substructures as the resolu- 
tion of the FMap reconstruction is made finer. These peaks 



are persistent across the FM^p realisations, and largely show 
minimal B-mode contamination. 

The identified structures were localised by comparing 
their detected positions across the 16 aperture mass recon- 
structions, and the 2a dispersion in these measurements was 
consistently less than ^ 10^^ in all cases, which highlights 
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the power of the flexion aperture mass statistic to precisely 
locate substructures in a cluster environment. 

By comparison with template models, a pseudo- 
probability distribution was computed for each peak in order 
to determine which mass models best fit the data. In each 
case, the distribution of allowed masses was fairly narrow, 
and all peaks showed consistency between the masses ob- 
tained using NFW and SIS template models. This implies 
that our technique offers a method to determine the masses 
of substructures using only the aperture mass reconstruc- 
tions and requiring no complicated parametric modelling. 

The results presented here show that Abell 1689 ex- 
hibits a rather complex internal structure. The presence of 
two dominant mass peaks (peaks [1,3] and 2) favour the idea 
that Abell 1689 has undergone a merger in its recent history 
(Andersson & Madejski, 2004); these two mass concentra- 
tions may be the remnants of this interaction. Our results 
also highlight that standard weak lensing methods might not 
provide a complete description of the structure of a cluster 
of galaxies; the shear aperture mass reconstructions appear 
to favour a single lensing halo, whilst the FM^p reconstruc- 
tions clearly show in much more detail the sub-arcminute 
structure of the cluster. 

Our method is seen to show excellent agreement with 
strong lensing mass models of the cluster near the central, 
critical region of the cluster. Moreover, FM^p techniques of- 
fer the possibility to detect significant substructure far re- 
moved from the critical region of the cluster, and where 
strong lensing measurements are not available. 

We have demonstrated that even in the absence of high 
signal to noise FM^p measurements, mass concentrations 
can be detected, and the allowed parameter space well con- 
strained, by using a moderate number of FM^p reconstruc- 
tions, and incorporating a requirement of persistence of the 
structures across reconstructions and no B-mode signal in 
order to avoid contamination from peaks arising out of noise. 
This technique offers a simple and straightforward method 
for mapping and quantifying the substructure content of 
clusters of galaxies without the need for dynamical infor- 
mation or parametric modelling. It will prove very valuable 
in the analysis of data such as the HST Multi-Cycle Trea- 
sury Program "CLASH" (P.I.: M. Postman and H. Ford), 
which will target 25 massive clusters. In combination with 
traditional weak lensing and strong lensing, flexion provides 
us with a multi-scale view of the distribution of dark matter 
on cluster scales. 
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Al Singular Isothermal Sphere Model 

We begin by considering the SIS model, in which the density 
profile is defined by: 



p{r) 



(Al) 



where ay is the three dimensional velocity dispersion of the 
lens. The convergence for this lens is given by: 



ni0) = 



(A2) 



where we have now made the transformation to angular co- 
ordinates = DdO), and the Einstein radius of the lens, 
Oe, is given by: 



Oe = 47r 



The deflection angle, shear, and flexions are 



a{0) 
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where is the position angle of the source with respect to 
the lens; i.e. = arctan(^2/^i). 



A2 Navarro-Frenk- White Profile 

This profile has a density defined as 

Scpc 



p{r) 



(r/r.)(l+r/r.)2 ' 
where pc is the critical density of the universe, and 
200 



6c 



3 ln(l+c) -c/(l + c) 



(A8) 



(A9) 



c is the concentration parameter, defined as i^2oo/Ts. Before 
defining the convergence, it is useful to introduce some short- 
hand notation. First, we define x = \^\/rs = |^|/^s, where 
1^1 = l^l/Dd- Second, we define a normalisation factor: 



2pc6crs 



(AlO) 



APPENDIX A: CONVERGENCE, SHEAR AND 
FLEXION PROFILES 

This appendix presents the expected shear and flexion pro- 
files arising from two models commonly used in gravitational 
lensing: the singular isothermal sphere (SIS) model and the 
Navarro-Frenk- White (NFW) profile. Lasky & Fluke (2009) 
provide a derivation of the expected shear and flexion signals 
from both of these profiles, as well as from a Sersic profile, 
and the results presented here are based in part on their 
work. 



In this notation, the convergence is given by 



where 



Eix) = { 



k{x) — 
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[1 - S(x) 
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The lensing operators are 



a{x) = 

7(x) = 

^{x) = 

Q{x) = 
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l-E(.)-2(l-^) (ln(|)+H(.)) 
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8 (^x - i j In (I) + 3(1 - 2x^) + (^15x^ -20 + -^^ E{x) 
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(A14) 
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(A16) 



where, again, is the position angle of the source with re- 
spect to the centre of the lens. Note that in the limit of 
X — 1, the lensing operators become: 
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